Downloading and cropping EJScreen data

Load the required R packages (install first as needed, tbeptools installation instructions are here) .

library(tbeptools)
library(sf)

Linking to GEOS 3.10.2, GDAL 3.4.1, PROJ 8.2.1; sf_use_s2() is TRUE

library(mapview)
library(dplyr)

Attaching package: 'dplyr'

The following objects are masked from 'package:stats':

    filter, lag

The following objects are masked from 'package:base':

    intersect, setdiff, setequal, union

library(RColorBrewer)

Download the relevant file from EJScreen. All files available at https://gaftp.epa.gov/EJSCREEN/2022/?C=S;O=A. The file is downloaded to a temporary directory.

# url with zip gdb to download
urlin <- 'https://gaftp.epa.gov/EJSCREEN/2022/EJSCREEN_2022_Supplemental_with_AS_CNMI_GU_VI_Tracts.gdb.zip'

# download file
tmp1 <- tempfile(fileext = "zip")
download.file(url = urlin, destfile = tmp1)

Unzip the geodatabase that was downloaded to a second temporary directory.

# unzip file
tmp2 <- tempfile()
utils::unzip(tmp1, exdir = tmp2, overwrite = TRUE)

Read the polygon layer from the geodatabase.

# get the layers from the gdb
gdbpth <- list.files(tmp2, pattern = '\\.gdb$', full.names = T)
lyrs <- st_layers(gdbpth)$name

# read the layer
dat <- st_read(dsn = gdbpth, lyrs)

Reading layer `EJSCREEN_Full_with_AS_CNMI_GU_VI' from data source 
  `/tmp/Rtmpu9wVd2/fileb64641422416/EJSCREEN_2022_Supplemental_with_AS_CNMI_GU_VI_Tracts.gdb' 
  using driver `OpenFileGDB'
Simple feature collection with 86000 features and 158 fields
Geometry type: MULTIPOLYGON
Dimension:     XY
Bounding box:  xmin: -19951910 ymin: -1617130 xmax: 16259860 ymax: 11554350
Projected CRS: WGS 84 / Pseudo-Mercator

Intersect the layer with the Tampa Bay watershed.

# intersect the layer with the tb watershed
dattb <- dat %>% 
  st_transform(crs = st_crs(tbshed)) %>% 
  st_make_valid() %>% 
  st_intersection(tbshed)

old-style crs object detected; please recreate object with a recent sf::st_crs()
old-style crs object detected; please recreate object with a recent sf::st_crs()

Warning: attribute variables are assumed to be spatially constant throughout
all geometries

View the data using mapview.

# view the layer
mapview(dattb)

dattbcopy <- dattb %>%
  filter(ACSTOTPOP > 0) %>%
  mutate(threshold_income = ifelse(P_LWINCPCT >= 80, 1, 0),
         threshold_unempl = ifelse(P_UNEMPPCT >= 80, 1, 0),
         threshold_lingui = ifelse(P_LNGISPCT >= 80, 1, 0),
         threshold_educat = ifelse(P_LESHSPCT >= 80, 1, 0),
         threshold_lifexp = ifelse(P_LIFEEXPCT >= 80, 1, 0)) %>%
  rowwise() %>%
  mutate(threshold_N = sum(threshold_income,threshold_unempl,threshold_lingui,threshold_educat,threshold_lifexp, na.rm = TRUE)) %>%
  mutate(underserved = ifelse(threshold_N > 1, "Yes", "No"))
mapview(dattbcopy, zcol = "threshold_N")

mapview(dattbcopy, zcol = "underserved")

The layer can be saved as an RData object if needed. Size is very small (~1mb).

# save the layer as an RData object (~1mb)
save(dattb, file = 'data/dattb.RData')

Unlink the temporary files so they’re deleted when you’re done.

# remove temp files
unlink(tmp1, recursive = TRUE)
unlink(tmp2, recursive = TRUE)